/*
  pyivl.c  -  Python bindings to IVL Library

  This file is part of the IVL Library (LIBIVL).
  
  LIBIVL - library for interval computations 
  
  Copyright (c) 2007 Mateusz Paprocki
  
  The IVL Library is free software; you can redistribute it and/or modify
  it under the terms of the GNU Lesser General Public License as published
  by the Free Software Foundation; either version 2.1 of the License, or
  (at your option) any later version.

  The IVL Library is distributed in the hope that it will be useful, but
  WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
  or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
  License for more details.

  You should have received a copy of the GNU Lesser General Public License
  along with the IVL Library; see the file LICENCE.  If not, write to
  the Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,
  Boston, MA 02110-1301, USA.
*/

#include <Python.h>

#include "interval.h"
#include "pyivl.h"

static char interval_doc [ ] = "\n";

static char pyivl_doc [ ] = "\n";




static char make_doc [ ] = 
	"make (lower, upper) -> interval\n";
	
static PyObject *interval_make(PyIntervalObject *self, PyObject *args) {
	
}

static char expand_doc [ ] =
	"expand (point, domain) -> interval\n";
	
static PyObject *interval_expand(PyIntervalObject *self, PyObject *args) {
	/*
	PyObject *object = NULL;
	int result;
	
	if (PyArg_ParseTuple(args, "O", &object))
		PyErr_SetString(PyExc_ValueError, "invalid argument");
	else {
		if (PyIvl_Check(object) 		
			result = ivl_strictly_inside(PyIvl_AsIVL(self), PyIvl_AsIVL(object));
		else if (PyInt_Check(object))
			result = ivl_strictly_inside_s(PyIvl_AsIVL(self), PyFloat_AS_LONG(object));
		else if (PyFloat_Check(object))
			result = ivl_strictly_inside_d(PyIvl_AsIVL(self), PyFloat_AS_DOUBLE(object));
#ifndef PYIVL_USE_GMPY
		else if (Pympz_Check(object))
			result = ivl_strictly_inside_z(PyIvl_AsIVL(self), Pympf_AS_MPZ(object));
		else if (Pympq_Check(object))
			result = ivl_strictly_inside_q(PyIvl_AsIVL(self), Pympf_AS_MPQ(object));
		else if (Pympf_Check(object))
			result = ivl_strictly_inside_f(PyIvl_AsIVL(self), Pympf_AS_MPF(object));
#endif /* not PYIVL_USE_GMPY */
		else {
			PyErr_SetString(PyExc_ValueError, "unsupported type");
			return NULL;
		}
		
		Py_DECREF(object);
		
		if (result) {
			Py_INCREF(Py_True);
			object = Py_True;
		} else {
			Py_INCREF(Py_False);
			object = Py_False;
		}
	}

	return object;
	*/
}

static char has_zero_doc [ ] =
	"has_zero () -> boolean\n\n"
	"Return True if zero is included within lower and upper bounds.\n"
	"In special case one or both of the bounds can be zero and this\n"
	"function will still return True. This is a shorter version of\n"
	"is_inside(interval, 0) invocation. Otherwise or when the domain\n"
	"is undefined return False.\n";
	
static PyObject *interval_has_zero(PyIntervalObject *self, PyObject *args) {
	if (ivl_has_zero(PyIvl_AsIVL(self))) {
		Py_INCREF(Py_True);
		return Py_True;
	} else {
		Py_INCREF(Py_False);
		return Py_False;
	}
}

static char is_strictly_inside_doc [ ] =
	"is_strictly_inside (point) -> boolean\n\n"
	"Return True if specified point exactly greater than the lower\n"
	"bound and exactly less than the upper bound. Otherwise or when\n"
	"the domain is undefined return False.\n";

static PyObject *interval_is_strictly_inside(PyIntervalObject *self, PyObject *args) {
	PyObject *object = NULL;
	int result;
	
	if (PyArg_ParseTuple(args, "O", &object))
		PyErr_SetString(PyExc_ValueError, "invalid argument");
	else {
		if (PyIvl_Check(object) 		
			result = ivl_strictly_inside(PyIvl_AsIVL(self), PyIvl_AsIVL(object));
		else if (PyInt_Check(object))
			result = ivl_strictly_inside_s(PyIvl_AsIVL(self), PyFloat_AS_LONG(object));
		else if (PyFloat_Check(object))
			result = ivl_strictly_inside_d(PyIvl_AsIVL(self), PyFloat_AS_DOUBLE(object));
#ifndef PYIVL_USE_GMPY
		else if (Pympz_Check(object))
			result = ivl_strictly_inside_z(PyIvl_AsIVL(self), Pympf_AS_MPZ(object));
		else if (Pympq_Check(object))
			result = ivl_strictly_inside_q(PyIvl_AsIVL(self), Pympf_AS_MPQ(object));
		else if (Pympf_Check(object))
			result = ivl_strictly_inside_f(PyIvl_AsIVL(self), Pympf_AS_MPF(object));
#endif /* not PYIVL_USE_GMPY */
		else {
			PyErr_SetString(PyExc_ValueError, "unsupported type");
			return NULL;
		}
		
		Py_DECREF(object);
		
		if (result) {
			Py_INCREF(Py_True);
			object = Py_True;
		} else {
			Py_INCREF(Py_False);
			object = Py_False;
		}
	}

	return object;
}

static char is_inside_doc [ ] =
	"is_inside (point) -> boolean\n\n"
	"Return True if specified point is included within lower and upper\n"
	"bounds. The point can be equal to one or both of the bounds. \n"
	"Otherwise or when the domain is undefined return False.\n";
	
static PyObject *interval_is_inside(PyIntervalObject *self, PyObject *args) {
	PyObject *object = NULL;
	int result;
	
	if (PyArg_ParseTuple(args, "O", &object))
		PyErr_SetString(PyExc_ValueError, "invalid argument");
	else {
		if (PyIvl_Check(object) 		
			result = ivl_inside(PyIvl_AsIVL(self), PyIvl_AsIVL(object));
		else if (PyInt_Check(object))
			result = ivl_inside_s(PyIvl_AsIVL(self), PyFloat_AS_LONG(object));
		else if (PyFloat_Check(object))
			result = ivl_inside_d(PyIvl_AsIVL(self), PyFloat_AS_DOUBLE(object));
#ifndef PYIVL_USE_GMPY
		else if (Pympz_Check(object))
			result = ivl_inside_z(PyIvl_AsIVL(self), Pympf_AS_MPZ(object));
		else if (Pympq_Check(object))
			result = ivl_inside_q(PyIvl_AsIVL(self), Pympf_AS_MPQ(object));
		else if (Pympf_Check(object))
			result = ivl_inside_f(PyIvl_AsIVL(self), Pympf_AS_MPF(object));
#endif /* not PYIVL_USE_GMPY */
		else {
			PyErr_SetString(PyExc_ValueError, "unsupported type");
			return NULL;
		}
		
		Py_DECREF(object);
		
		if (result) {
			Py_INCREF(Py_True);
			object = Py_True;
		} else {
			Py_INCREF(Py_False);
			object = Py_False;
		}
	}

	return object;
}

static char is_point_doc [ ] =
	"is_point () -> boolean\n\n"
	"Return True if lower and upper bounds are equal or are\n"
	"infinities of the same sign. Otherwise or when the domain\n"
	"is undefined return False.\n";
	
static PyObject *interval_is_point(PyIntervalObject *self, PyObject *args) {
	if (ivl_point_p(PyIvl_AsIVL(self))) {
		Py_INCREF(Py_True);
		return Py_True;
	} else {
		Py_INCREF(Py_False);
		return Py_False;
	}
}

static char is_infinity_doc [ ] =
	"is_infinity () -> boolean\n\n"
	"Return True if lower and upper bounds are infinities of the same\n"
	"sign. Otherwise or when the domain is undefined return False.\n";
	
static PyObject *interval_is_infinity(PyIntervalObject *self, PyObject *args) {
	if (ivl_point_inf_p(PyIvl_AsIVL(self))) {
		Py_INCREF(Py_True);
		return Py_True;
	} else {
		Py_INCREF(Py_False);
		return Py_False;
	}
}

static char is_bounded_doc [ ] =
	"is_bounded () -> boolean\n\n"
	"Return True if both lower and upper bounds are numerical values\n"
	"(not infinites), which don't have to be exactly the same. Otherwise\n"
	"or when the domain is undefined return False.\n";
	
static PyObject *interval_is_bounded(PyIntervalObject *self, PyObject *args) {
	if (ivl_bounded_p(PyIvl_AsIVL(self))) {
		Py_INCREF(Py_True);
		return Py_True;
	} else {
		Py_INCREF(Py_False);
		return Py_False;
	}
}

static char is_unbounded_doc [ ] =
	"is_unbounded () -> boolean\n\n"
	"Return True if lower and upper bounds are distinct infinities.\n"
	"Otherwise or when the domain is undefined return False.\n";
	
static PyObject *interval_is_unbounded(PyIntervalObject *self, PyObject *args) {
	if (ivl_unbounded_p(PyIvl_AsIVL(self))) {
		Py_INCREF(Py_True);
		return Py_True;
	} else {
		Py_INCREF(Py_False);
		return Py_False;
	}
}

/**
 *
 */
static PyMethodDef interval_methods [ ] = {
	{ "make", (PyCFunction)interval_make, METH_VARARGS, make_doc },
	{ "expand", (PyCFunction)interval_expand, METH_VARARGS, expand_doc },

	{ "has_zero", (PyCFunction)interval_has_zero, METH_NOARGS, has_zero_doc },
	
	{ "is_strictly_inside", (PyCFunction)interval_is_strictly_inside, METH_VARARGS, is_strictly_inside_doc },
	{ "is_inside", (PyCFunction)interval_is_inside, METH_VARARGS, is_inside_doc },
	
	{ "is_point", (PyCFunction)interval_is_point, METH_NOARGS, is_point_doc },
	{ "is_infinity", (PyCFunction)interval_is_infinity, METH_NOARGS, is_infinity_doc },
	{ "is_bounded", (PyCFunction)interval_is_bounded, METH_NOARGS, is_bounded_doc },
	{ "is_unbounded", (PyCFunction)interval_is_unbounded, METH_NOARGS, is_unbounded_doc },
	
	{ NULL, NULL, 0, NULL }		/* sentinel */
};

// getsetters : prec, domain
// as_seq : length -> len

/** Make new interval object. This function will allocate object's memory
 *  and initialize interval data, which means lowe un upper bunds are set
 *  to NaN and the domain is set to IVL_UNDEFINED. The default precision
 *  is used for new interval. This can be overriden by setting explicitely
 *  new precision using ivl_set_default_prec module function. In fact,
 *  setting precision at this moment is quite useless, so one can use
 *  __init__ function and pass the desired precision in it. 
 */
static PyObject *interval_new(PyTypeObject *type, PyObject *args, PyObject *kwds) {
    PyObject *self;

    self = type->tp_alloc(type, 0);
	
    ivl_init(((PyIntervalObject *)self)->ob_ivl); /* use default precision */

    return (PyObject*)self;
}

/** Dispach given Python object through possible types to interval functions,
 *  which include integer, double, string, mpz, mpq and mpf. If one of the
 *  assertions is successful then coerced value in MPFR format is returned,
 *  rounded in specified direction. Otherwise an arbitrary expression is
 *  executed, which will possibly goto somewhere or rise an exception.
 */
#ifndef PYIVL_USE_GMPY
	#define PYIVL_DISPATCH_ARG(object, fr, rnd, do_else)           \
		if (PyInt_Check(object)) {                                 \
			mpfr_set_si(fr, PyInt_AS_LONG(object), rnd);           \
		} else if (PyFloat_Check(object)) {                        \
			mpfr_set_d(fr, PyFloat_AS_DOUBLE(object), rnd);        \
		} else if (PyString_Check(object)) {                       \
			mpfr_set_str(fr, PyString_AS_STRING(object), 10, rnd); \
		} else {                                                   \
			do_else;                                               \
		}
#else
	#define PYIVL_DISPATCH_ARG(object, fr, rnd, do_else)           \
		if (PyInt_Check(object)) {                                 \
			mpfr_set_si(fr, PyInt_AS_LONG(object), rnd);           \
		} else if (PyFloat_Check(object)) {                        \
			mpfr_set_d(fr, PyFloat_AS_DOUBLE(object), rnd);        \
		} else if (PyString_Check(object)) {                       \
			mpfr_set_str(fr, PyString_AS_STRING(object), 10, rnd); \
		} else if (Pympz_Check(object)) {                          \
			mpfr_set_z(fr, Pympz_AS_MPZ(object), rnd);             \
		} else if (Pympq_Check(object)) {                          \
			mpfr_set_q(fr, Pympq_AS_MPQ(object), rnd);             \
		} else if (Pympf_Check(object)) {                          \
			mpfr_set_f(fr, Pympf_AS_MPF(object), rnd);             \
		} else {                                                   \
			do_else;                                               \
		}
#endif /* not PYIVL_USE_GMPY */

/** Initialize interval with new data. There are to schemes. User can specify
 *  an existing interval, so its bounds and domain will be copied. Or both
 *  lower and upper bounds can be specified explicitely as __init__ arguments.
 *  Supported are all types accepted by IVL Library (see PYIVL_DISPATCH_ARG).
 *  An additional argument can be specified to set desired precision. If not
 *  specified input interval precision will be used or MPFR defaults.
 */
static int interval_init(PyIntervalObject *self, PyObject *args, PyObject *kwds) {
	static char *kwlist [ ] = { "precision", NULL };
	PyObject *o_lower, *o_upper, *ivl;
	mpfr_t f_lower, f_upper;
	mp_prec_t prec = 0;

	if (PyArg_ParseTupleAndKeywords(args, kwds, "O|i", kwlist, &ivl, &prec)) {
		if (PyInterval_Check(ivl)) {
			if (prec != 0) 
				ivl_set_prec(self->ob_ivl, prec);
			else
				ivl_set_prec(self->ob_ivl, ivl_get_prec(((PyIntervalObject *)ivl)->ob_ivl));
			
			ivl_set(self->ob_ivl, ((PyIntervalObject *)ivl)->ob_ivl);			
		} else {
			if (prec != 0) 
				ivl_set_prec(self->ob_ivl, prec);
			
			if (PyInt_Check(ivl)) {
				ivl_set_s(self->ob_ivl, PyInt_AS_LONG(ivl));
			} else if (PyFloat_Check(o_lower)) {
				ivl_set_d(self->ob_ivl, PyFloat_AS_DOUBLE(ivl));
			} else if (PyString_Check(o_lower)) {
				ivl_set_str(self->ob_ivl, PyString_AS_STRING(ivl), 10);
#ifdef PYIVL_USE_GMPY
			} else if (Pympz_Check(o_lower)) {
				ivl_set_z(self->ob_ivl, Pympz_AS_MPZ(ivl));
			} else if (Pympq_Check(o_lower)) {
				ivl_set_q(self->ob_ivl, Pympq_AS_MPQ(ivl));
			} else if (Pympf_Check(o_lower)) {
				ivl_set_f(self->ob_ivl, Pympf_AS_MPF(ivl));
#endif /* PYIVL_USE_GMPY */
			} else {
				Py_DECREF(ivl);
				goto bounds;
			}
			
			if (!ivl_set_domain(self->ob_ivl, IVL_TOTAL)) {
				PyErr_SetString(PyExc_ValueError, "one of the bounds is NaN");
				
				Py_DECREF(ivl);
				return -1;
			}
		}
		
		Py_DECREF(ivl);
		return 0;
	}
	
bounds:	
	if (!PyArg_ParseTupleAndKeywords(args, kwds, "OO|i", kwlist, &o_lower, &o_upper, &prec))
		return -1;

	if (prec == 0) {
		mpfr_init(f_lower); mpfr_init(f_upper);
	} else {
		ivl_set_prec(self->ob_ivl, prec);
		mpfr_init2(f_lower, prec);
		mpfr_init2(f_upper, prec);
	}

	PYIVL_DISPATCH_ARG(o_lower, f_lower, GMP_RNDD, goto error)
	PYIVL_DISPATCH_ARG(o_upper, f_upper, GMP_RNDU, goto error)

	Py_DECREF(o_lower);
	Py_DECREF(o_upper);

	ivl_make_f(self->ob_ivl, f_lower, f_upper);

	mpfr_clear(f_lower);
	mpfr_clear(f_upper);

	return 0;

error:
	PyErr_SetString(PyExc_ValueError, "bounds have invalid type");

	mpfr_clear(f_lower);
	mpfr_clear(f_upper);

	Py_DECREF(o_lower);
	Py_DECREF(o_upper);

	return -1;
}

/** Deallocate interval object. This function will only clear the interval, 
 *  free memory occupied by its data and finally free object's memory.  
 */
static void interval_dealloc(PyObject *self) {
	ivl_clear(((PyIntervalObject *)self)->ob_ivl);    
    self->ob_type->tp_free((PyObject*)self);
}


static int interval_print(PyIntervalObject *self, FILE *fp, int flags) {
	ivl_print(fp, self->ob_ivl, 10, 0);
	return 0;
}

static PyObject *interval_repr(PyIntervalObject *self) {
	char buf[100];
	
	return PyString_FromString(buf);
}

static PyObject *interval_str(PyIntervalObject *self) {
	char buf[100];
	
	return PyString_FromString(buf);
}

static PyObject *interval_richcompare(PyObject *x, PyObject *y, int op) {
	
}

/** Compute hash value for eq. dictionery use. For now the algorithm is to
 *  take hash value of both bounds coerced to doubles (with respect to proper
 *  inifinities) and to make exclusive-or on those values.
 *
 *  TODO : The algorithm used here is trivial and it does not take domain
 *  into account. There should be done future improvements. 
 */
static long interval_hash(PyIntervalObject *self) {
	long hash_left, hash_right;
	
	hash_left = _Py_HashDouble(ivl_get_left_d(self->ob_ivl));
	hash_right = _Py_HashDouble(ivl_get_right_d(self->ob_ivl));
	
	return hash_left ^ hash_right;
}

/** This defines interval as a Python base object, which means it is first 
 *  in inheritance scheme and, at least for now, it does not allow any
 *  sub-classing. Intervals implement number protocol including persistent
 *  and in-place operations. 
 */ 
PyTypeObject PyInterval_Type = {
    PyObject_HEAD_INIT(NULL)
    0,											/* ob_size */
    "interval",        							/* tp_name */
    sizeof(PyIntervalObject),					/* tp_basicsize */
    0,											/* tp_itemsize */
    (destructor)interval_dealloc,				/* tp_dealloc */
    (printfunc)interval_print,					/* tp_print */
    0,											/* tp_getattr */
    0,											/* tp_setattr */
    0,											/* tp_compare */
    (reprfunc)interval_repr,					/* tp_repr */
	0,
    //&interval_as_number,						/* tp_as_number */
    0,											/* tp_as_sequence */
    0,											/* tp_as_mapping */
    (hashfunc)interval_hash,					/* tp_hash */
    0,											/* tp_call */
    (reprfunc)interval_str,						/* tp_str */
    0,											/* tp_getattro */
    0,											/* tp_setattro */
    0,											/* tp_as_buffer */
    Py_TPFLAGS_DEFAULT,						 	/* tp_flags */
    interval_doc,								/* tp_doc */
    0,											/* tp_traverse */
    0,											/* tp_clear */
    interval_richcompare,						/* tp_richcompare */
    0,											/* tp_weaklistoffset */
    0,											/* tp_iter */
    0,											/* tp_iternext */
    interval_methods,							/* tp_methods */
    0,											/* tp_members */
    0,											/* tp_getset */
    0,											/* tp_base */
    0,											/* tp_dict */
    0,											/* tp_descr_get */
    0,											/* tp_descr_set */
    0,											/* tp_dictoffset */
    (initproc)interval_init,					/* tp_init */
    0,											/* tp_alloc */
    interval_new,								/* tp_new */
};

/**
 *
 */
static PyMethodDef pyivl_methods [ ] = {
	{ NULL, NULL, 0, NULL }		/* sentinel */
};

/** Module initialization function
 *
 *  This will initialize PyIVL module and add new objects to it. For now
 *  there is only one - PyIntervalObject, but in the future complex and
 *  multidimensional intervals will be added.
 *
 *  If specified it will also import Gmpy library, which will make types
 *  like mpz, mpq and mpf for unbounded computations available. 
 */
PyMODINIT_FUNC initpyivl(void) {
    PyObject *module = NULL;

    if (PyType_Ready(&PyInterval_Type) < 0)
        return;

    if ((module = Py_InitModule3("pyivl", pyivl_methods, pyivl_doc)) == NULL)
		return;

    Py_INCREF(&PyInterval_Type);
		
    PyModule_AddObject(module, "interval", (PyObject *)&PyInterval_Type);
}
